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. Abstract We consider non-reversible perturbations of reversible diffusions 

■ that do not alter the invariant distribution and we ask whether there exists 

an optimal perturbation such that the rate of convergence to equilibrium is 
maximized. We solve this problem for the case of linear drift by proving the 
existence of such optimal perturbations and by providing an easily imple- 
i mcntable algorithm for constructing them. We discuss in particular the role 

of the prefactor in the exponential convergence estimate. Our rigorous results 
are illustrated by numerical experiments. 
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1 Introduction 



> 

■ 1.1 Motivation 

p: 

. The problem of convergence to equilibrium for diffusion processes has at- 

tracted considerable attention in recent years. In addition to the relevance of 
P\J . this problem for the convergence to equilibrium of some systems in statistical 
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physics, see for example [29], such questions are also important in statistics, 
for example in the analysis of Markov Chain Monte Carlo (MCMC) algo- 
rithms [S] . Roughly speaking, one measure of efficiency of an MCMC algo- 
rithm is its rate of convergence to equilibrium, and increasing this rate is 
thus the aim of many numerical techniques (see for example [5]). 

Let us recall the basic approach for a reversible diffusion. Suppose that 
we are interested in sampling from a probability distribution function 

^ f e-y dx ' 

where V : — >■ M is a given smooth potential such that J^j^ e^^ dx < oo . 
A natural dynamics to use is the reversible dynamics 

dXt^ -VViXt)dt + V2dWt, (2) 

where Wt denotes a standard A^-dimensional Brownian motion. Let us denote 
by ipt the probability density function of the process Xt at time t. It satisfies 
the Fokker-Planck equation 

dtiJt^V -ivv^t+y^t) ■ (3) 

Under appropriate assumptions on the potential V (e.g. that i|VT^(x)p — 
AV{x) — > -l-oo as — ^ -l-oo , see dU A. 19]), the density -(/'oo satisfies a 
Poincare inequality: there exists A > such that for all probability density 
functions 6, 
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1 ipoodx < 

Woo J A 



2 

iJoodx . (4) 



The optimal parameter A in is the opposite of the smallest (in absolute 
value) non-zero eigenvalue of the Fokker-Planck operator V • {W ■ -I-V-), 
which is self-adjoint in L^(K^, dx) (see ([7]) below). Thus, A is also called 
the spectral gap of the Fokker-Planck operator. 

It is then standard to show that (H]) is equivalent to the following inequal- 
ity, which shows exponential convergence to the equilibrium for ([2|): for all 
initial conditions -00 G ^i^i^-^ I'^'Zo'dx), for all times i > 0, 

llV't - V'oo|li2(^-l) < e^^*||V0 - V'oo|Il2(^-1) , (5) 

where || • 11^2(^-1) denotes the norm in (M^ , ■0- 1 ) , namely \\f\\\^^^_,^ = 

Jt^n P{x)il'^{x) dx . This equivalence is a simple consequence of the follow- 
ing identity: if tpt is solution to then 

d /■ / ; \ 2 

^ll^t-V'oo|!i2(^-l)--2 



■0OO 



ipoodx . (6) 



In view of ([5]), the algorithm is efficient if A is large, which is typically not 
the case if Xt is a metastable process (see [25]). A natural question is there- 
fore how to design a Markovian dynamics which converges to equilibrium 
distribution ipao (much) faster than ^ . There are many approaches (impor- 
tance sampling methods, constraining techniques, see for example |26j). and 
the focus here is on modifying the dynamics ([2]) to a non-reversible dynamics, 
which has the same invariant measure. 
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1.2 Non- reversible diffusion 

As noticed in p51l24| . one way to accelerate the convergence to equilibrium 
is to depart from reversible dynamics (see also |10| for related discussions 
for Markov Chains). Let us recall that the dynamics ([2]) is reversible in the 
sense that if Xq is distributed according to ipoo(x)dx, then {Xt)o<t<T and 
{XT-t)Q<t<T have the same law. This is equivalent to the fact that the 
Fokker-Planck operator is self-adjoint in L'^{M.^ ,ip^^dx): 

V • (VFV + VV) Voo'c^a: = - / V {H^^') ■ V {^^-') i^oodx 

V • (yV(i3 + V(i3)ip^;^^dx. 

Now, a natural non-reversible dynamics to sample from the distribution 
4,^{x)dx= dx is: 

dX^^ = ( - "^ViX^) + b{X^)) dt + V2dWt, (8) 

where b is taken to be divergence-free with respect to the invariant distribu- 
tion il^aoix) dx: 

V • (6e-^) = 0, (9) 

so that "000 {x) dx is still the invariant measure of the dynamics ([5]) . A general 
way to construct such a 6 is to consider 

b = JVV, (10) 

where J is a constant antisymmetric matrix. 

It is important to note that the dynamics (O is non-reversible. Indeed, 
one can check that (X('')o<t<T has the same law as {X^^^)f)^t<T (notice 
the minus sign in front of 6), and thus not the same law as (Ary_j)o<t<T- 
Likewise, Equation ([7]) now becomes: 

V-((Vy-6)V + V0)0i/;-Ma;= f V ■ {{VV + b)(l) + Vc/)) ^ i^-^dx . 

Again, notice the change of sign in front of b. 

From ITU]) it is clear that there are many (in fact, infinitely many) differ- 
ent ways to modify the reversible dynamics without changing the invariant 
measure. A natural question is whether the addition of a non-reversible term 
can improve the rate of convergence to equilibrium and, if so, whether there 
exists an optimal choice for the perturbation that maximizes the rate of 
convergence to equilibrium. The goal of this paper is to present a complete 
solution to this problem when the drift term in ^ is linear. 

More precisely, let denote the law of the process X^, i.e. the solution 
to the Fokker-Planck equation 



(12) 
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Using the fact that ipoo is a stationary solution to (fT2)) (which is equivalent 
to (O) and under the assumption that tpoo satisfies the Poincare inequal- 
ity (|H), one can check that the upper bound for the reversible dynamics ^ 
is still valid: 

y4 - V'oolL2(^^i) < e-^*||Vg - V'oolL2(^-i) . (13) 

Actually, as in the reversible case, (IT^ (for all initial conditions ipo) is equiv- 
alent to dH). This is because ([6]) also holds for ip'' solution to (fT2)) . In other 
words, adding a non-reversible part to the dynamics cannot be worse than 
the original dynamics ^ (where & = 0) in terms of exponential rate of con- 
vergence. 

What we show below (for a linear drift) is that it is possible to choose b 
in order to obtain a convergence at exponential rate of the form: 

Ut - V^OO 1^2(^-1) < C{V, b)e-^'\\4 - V-OO 1^2(^-1) , (14) 

with A > A and C{V, b) > 1. It is important to note the presence of the 
constant C{V,b) in the right-hand side of ((T4l) . For a reversible diffusion 
(6 = 0), the spectral theorem forces the optimal C{V,0) to be equal to one, 
and A = A, the Poincare inequality constant of V'oo (since ([5]) implies (H])). 
The interest in adding a non-reversible perturbation is precisely to allow for 
a constant C{V,b) > 1 , which permits a rate A > A . The difficulty is thus to 
design a b such that A is large and C{V, b) is not too large. In the following, 
we adapt a two-stage strategy: we first optimize b in order to get the largest 
possible A, and then we discuss how the constant C(V, b) behaves for this 
optimal rate of convergence. 

1.3 Bibliography 

This problem was studied in [23j for a linear drift (namely V is quadratic and 
b is linear) and in [21] for the general case. It was shown in these works that 
the addition of a drift function b satisfying © helps to speed up convergence 
to equilibrium. Furthermore, the optimal convergence rate was obtained for 
the linear problem (see also Proposition [T] in the present paper) and some 
explicit examples were presented, for ordinary differential equations in two 
and three dimensions. 

The behavior of the generator of the dynamics under a strong non- 
reversible drift has also been studied [4ll6lll3j. It was shown in [13 that the 
spectral gap attains a finite value in the limit as the strength of the perturba- 
tion becomes infinite if and only if the operator b ■ V has no eigenfunctions in 
an appropriate Sobolev space of index 1. These works, although relevant to 
our work, are not directly related to the present paper since our main focus is 
in obtaining the optimal perturbation rather than an asymptotic result. The 
effect of non-reversible perturbations to the constant in logarithmic Sobolev 
inequalities (LSI) for diffusions have also been studied, see In this paper 
examples were presented where the addition of a non-reversible perturbation 
can improve the constant in the LSI. 
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This work is also related to [15], where the authors use another idea to 
enhance the convergence to equilibrium. The principle is to keep a reversible 
diffusion, but to change the underlying Riemannian metric by considering 

dXf^' = -DW{X^') dt + VwdWt 

for a well chosen matrix D. More precisely, the authors apply this technique 
to a Hybrid Monte Carlo scheme. It would be interesting to set up some test 
cases in order to compare the two approaches: non-reversible drift versus 
change of the underlying metric. 

Finally, we would like to mention related recent works on spectral prop- 
erties of non-selfadjoint operators see for example [TI IHIHU] and references 
therein. 



1.4 Outline of the paper 

In this paper, we study the case of a linear drift. Namely, we consider ([2]) 
with a quadratic potential 

V{x) = '^x'^Sx, (15) 

where S* is a positive definite N x N symmetric matrix. In the following, 
we denote 5Ar(M) the set of symmetric matrices and 5^'^(M) the set of pos- 
itive definite symmetric matrices. The equilibrium distribution thus has the 
density 

, , , det(5)i/2 / x^Sx\ 

^^(^) - (2^)^/2 [ —) ■ (16) 

It can be checked that if the vector field b{x) is linear, it satisfies © if and 
only if 6 = —JSx with J — —J^ an antisymmetric real matrix, see Lemma [T] 
For a given S, the question is thus how to choose J in order to optimize the 
rate of convergence to equilibrium for the dynamics ([8]), which in our case 
becomes: 

dXt' ^ -{I + J)SXt' dt + V2dWt, (17) 

where / denotes the identity matrix in A4 n (IR) , the set oi N x N real valued 
matrices. 

We provide an answer to this question. In particular: 

1. We prove that it is possible to build an optimal J (denoted Jopt), which 
yields the best possible rate A (denoted Xopt) in ([T^. 

2. We provide an algorithm for constructing an optimal matrix Jopt- 

3. We obtain estimates on the constant C{V, h) = C{S, J) in ([Tl)) . 

It appears that this procedure becomes particularly relevant in the situation 
when the condition number of S is large (namely for an original dynamics 
with multiple timescales, see Sections 13.31 and [5)1. Discussions about the size 
of C(y,&) with respect to this conditioning and to the dimension N can be 
carried out very accurately. 
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The reason why the case of hnear drift is amenable to analysis is because 
it can be reduced to a linear algebraic problem, at least for the calculation of 
Xopt and the construction of Jopt ■ One way to understand this is the following 
remark: the spectrum of an operator of the form (which is precisely the form 
of the generator of the dynamics ((TT)) ) 

C = -{Bx) - V + A, (18) 

can be computed in terms of the eigenvalues of the matrix B. Here, B denotes 
any real square matrix with positive real spectrum. In [30 (see also 
and Proposition 1101 below) . it was indeed proven that the spectrum of C in 
LP spaces weighted by the invariant measure of the dynamics {p > 1) consists 
of integer linear combinations of eigenvalues of B: 

<y{^) = ^-pn,X,,n,cny (19) 

where {Ajj^^j^ denote the r (distinct) eigenvalues of B. In particular, the 
spectral gap of the generator C is determined by the eigenvalues of B, and 
this yields a simple way to design the optimal matrix Jopt- On the other 
hand, the control of the constant C{S, J) requires a more elaborate analysis, 
using Wick (in the sense of Wick ordered) calculus, see Section [Ql below. 

Compared to the related previous paper [33] , our contributions are three- 
fold: (i) we propose an algorithm to build the optimal matrix Jopt, (ii) we 
discuss how to get estimates on the constant C{S, J) and (iii) we consider 
the longtime behavior of the partial differential equation (|T2]) and not only 
of ordinary differential equations related to (jl7p . In particular, our analysis 
covers also non Gaussian initial conditions for the SDE (IT71) . Although the 
results that we obtain have a limited practical interest (there exist many 
efficient techniques to draw Gaussian random variables), we believe that this 
study is a first step towards further analysis, in particular for nonlinear drift 
terms. 

The rest of the paper is organized as follows. In Section [T75] we present the 
main results of this paper. In Section [2] we perform some preliminary calcu- 
lations. The linear algebraic problem and the evolution of the corresponding 
ordinary differential equation are studied in Section [31 Direct computations 
of the expectations and the variances are performed in Section|3]for Gaussian 
initial data. The convergence to equilibrium for the non-reversible diffusion 
process for general initial data is then studied in Section [S] Results of nu- 
merical simulations are presented in Section [51 Finally, some background 
material on Wick calculus, which is needed in the proofs of our main results, 
is presented in Appendix 1X1 

1.5 Main results 

For a potential given by (llSp . our first result is a simple lemma which char- 
acterizes all non-reversible perturbations that satisfy divergence-free condi- 
tion 
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Lemma 1 Let V{x) be given by (fTSj) and let b{x) = —Ax where A ^ M Ar(M). 
Then ^ is satisfied if and only if 

A = JS, with J ^-J^ . (20) 

Proof Equation ^ with b = — Ax and quadratic potential gives 



V 

which is equivalent to 



(^Axe-"^^ = 



Tr{A) + (AxfiSx) ^0 Vx G . 
This is equivalent to the conditions 

Tr(A) = and {A^ S) ^ -{A^ Sf . 
Set now J = AS^'^ . We have 

Tt{JS) = and S{J + J^)S = , 
which is equivalent to 

J=-J^. 

a 

We will denote the set of iV x TV real antisymmetric matrices by .4Ar(M) C 
Ai N (R) • The following result concerns the optimization of the spectrum of 
the matrix Bj = (/ + J)S, which appears in the drift of the dynamics P7)) 
and plays a crucial role in the analysis; see Equation (|19p . 

Theorem 1 Define Bj = {I + J)S . Then 

max minRe(cr(Bj)) = ^^4r^. (21) 

Furthermore, one can construct matrices J opt G -4jv(lK) such that the maxi- 
mum in (I2ip is attained. The matrix Jopt can be chosen so that the semi-group 
associated to Bj^ ^ satisfies the bound 



<&^'n{SY''e^v(-^t] , (22) 



where the matricial norm is induced by the euclidean norm on and k(S) 
\\S\\ ||S'^"'"|| denotes the condition number. 



Theorem [T] is a straightforward consequence of Proposition 0] and Proposi- 
tion [S] below, with an explicit expression for the constant C^'' given by ([50]). 

This expression allows to discuss the dependence of C^' on the dimension N 
(see Remark [5] for details). 

The partial differential equation version of this result requires to introduce 
the generator 

£j = -{Bjx)-\/ + A 
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7 ^oo^*^: 



of the semigroup {e ■')t>o considered in L 
(see (dll)), 

, , , det(5)i/2 / x'^Sx 

^-(-) = 75^<^-p [ — — 

Here L'^{M'^ ,ipaodx;C) denotes the set of functions / 
/rn \f\'^{x)ipoo{x)dx < oo. 



where, we recall 



pN 



C such that 



Theorem 2 For Bj = {I + J)S with J £ An , the drift- diffusion operator 
Cj = —{Bjx).W + A defined in L'^{M.^ ,ipoodx;C) with domain of definition 



D{Cj) = {we L^(M^,Voorfa;;C), Cju £ L^(R^ ,iP^dx;C)} 

generates a contraction semigroup {e*'^^')t>o o-nd it has a compact resolvent. 
Optimizing its spectrum with respect to J gives 



max minRc (ct(— £j) \ {0}) = 



Tr(^) 

N 



(23) 



Furthermore, the maximum in (j23p is attained for the matrices J opt G An ( 
constructed as in Theorem\^ The matrix J opt can be chosen so that 



uipoodx 



Tr(5) 
N 



u — 



UTpoodx 



(24) 



holds for all w G L^(l 
the condition number. 



'^'^ ,4'oodx;C) and all t > 0, where k{-) again denotes 



Theorem[5]is a straightforward consequence of PropositionlT^below. with an 
explicit expression for the constant C}^' given by Again, the dependence 

of C]y on the dimension N is discussed in Remark [TUl A simple corollary of 
this result is the following: 

Corollary 1 Let us consider the Fokker Planck equation associated to the 
dynamics (flT)) on X/; 



9tV/ = V-(i3ja;^/ + V^/) 



(25) 



where Bj = (/ + J)S. Let us assume that tp^ S L'^ (M.^ , tp^ dx) . Then, by 
considering J ~ — J opt, where J opt G -4Af(K) refers to the matrix considered 
in Theorem to get (j24p . Then the inequality 



Wt - V'oo 11^.(^-1) < CfniSf'^ exp 



Tr(g ) 

N 



holds for all t > , when V'oo *s defined by (fTB]l . 
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Proof This result is based on the following simple remark: i/)/ is a solu- 
tion to (EH) in L2(M^,-0^ida;) if and only if V/fc^ = e*^--'(-0^ -0^^) in 
L^(K^, ipao dx). Notice the minus sign in Then the exponential conver- 
gence is obtained from using the equality: 



□ 



Remark 1 A more general result but with a less accurate upper bound is 
given in Proposition [8l □ 

Remark 2 The partial differential equation 

dtf = Cjf = -{Bjx)-V! + Af 

which we consider in Theorem [5] is sometimes called the backward Kol- 
mogorov equation associated with the dynamics ([T7|. It is related to this 
stochastic differential equation through the Feynman-Kac formula: 

f{t,x)=Wf{X'^) 

where is the solution to (|17p and indicates that we consider a solution 
starting from x G M^: Xq = x. The partial differential equation 

which we consider in the Corollary [T] is the Fokker Planck (or forward Kol- 
mogorov) equation associated with (flT)) : if X^ ^ i'ii^) dx, then for all times 
t > 0, ipt is the probability density function of X^ . 

As explained in the proof of Corollary [1] above, these two partial differ- 
ential equations are related through a conjugation. See also, e.g. [11115^ . 

□ 

Remark 3 It would be interesting to explore extensions of this approach to 
the Langevin dynamics: 

J dqt = pt dt, 

\ dpt = -V\/(gt) dt - -ipt dt + y^dWt, 

which is ergodic with respect to the measure exp(— — |pp/2) dpdq. 
For example the following modification 

J dqt = (/ - J)pt dt, 

\ dpt - -(/ + J)VF(gt) dt - jpt dt + y^dWt, 

where J is an antisymmetric matrix leaves the measure exp{~-V{q) — 
\p\'^/2)dpdq stationary. In the linear case V{x) — ^ , this leads to a 
Kramers-Fokker-Planck operator which is a differential operator (at most) 
quadratic in {q,p,dq,dp) . Then the exponential decay rate can be reduced 
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to some (more involved) linear algebra problem following [20]. About the 
constant prefactor in front of the decaying in time exponential, the argu- 
ment based on sectiorality used in Lemma [3] does not apply anymore. It has 
to be replaced by hypoelliptic estimates in the spirit of pTifTOll^ . The ref- 
erence [20j provides accurate results for differential operators with at most 
quadratic symbols. □ 

Remark 4 We notice that the fundamental property div(6e^^) = is still 
satisfied for b{t,x) = J{t)VV{x), where J{t) is a time-dependent (determin- 
istic) antisymmetric matrix. This could be useful for further generalization 
of this approach. □ 



2 A useful rescaling 

The analysis will be carried out in a suitable system of coordinates which 
simplifies the calculations and the presentation of the intermediate results. 
We will perform one conjugation and a change of variables. 

First, from the partial differential equation point of view, it appears to 
be useful to work in L'^{m^,dx;C) instead of L'^{R^ ,ipoodx;C), since this 
allows to use standard techniques for the spectral analysis of partial differen- 
tial equations. In the following, the norm in L^(IR^, dx; C) is simply denoted 

II • 11^2 . For a general potential V, the mapping u i— > ^poo u sends unitar- 
ily L^(M^, dx; C) into L'^{M.^ ,ipoodx;C) with the associated transformation 
rules for the differential operators: 

V rj-\ V rj-t 1 T"! 

V^' = (a,,,...,a,„), (v^x-divx). 



is transformed into 

£ = e'^Ce^ = A- \\S/V\'^ + ]-AV + fo^V + ]:h^VV . (26) 
4 2 2 

In the linear case we consider in this paper, V{x) = ^x"^ Sx (where S — S"^ 
is positive definite), h{x) — —Ax and A = JS, J G ^Ar(]R), (see Lemma[T]), 
so that the operator 

C = Cj = -{BjxfV + A with Bj = {I + J)S, 




becomes 



Lj^A- -x'^S^x + -Tr(S') - x'^A^V - -x^A^Sx 
4 2 ^ ^ 2 

= A- -x'^S^x + -Tr(S') - x^SJ^S/ - -x'^Sj'^Sx 
4 2 ^ ^ 2 

= A- ^x'^S^x + ^Tr(S') -I- ^{x'^SJW - W^JSx) . 



11 



For the last line we have used 

= -J , x'^SJ^Sx = , 



0. 



with B^SJ^ , B^ ^ -JS and Tr(5'./) = Tr (^S^/^JS^/^^ 

According to Lemma [U we know that the kernel of Cj is Ce^^ — Ce^ " . 
The operator £j is unitarily equivalent to the operator Cj. 

The aim of the second change of variables is to modify the kernel of the 
operator to a centered Gaussian with covariance matrix being the identity. 
Let us introduce the new coordinates 



Then the operator Cj becomes: 



£j = V^^V, - ^y^Sy + ^TriS) + ^(j/^ JV, - Jy) 



(27) 



where 



The corresponding stochastic process is, in the new coordinate system {Yt = 

dYt = -{S + J)Ytdt + V2S^^^dWt . 

The L^-normalized element of ker >Cj is now simply the standard Gaussian 
distribution 

1 

(27r)^/4'' ' ■ 

Notice that Cj is still acting in L^(IR^, dx; C) . 
As a summary, u{t, x) satisfies 

dtu = Cju 

if and only if v(t, y) — ^/ipooiS^^^'^y) u{t, S^^^'^y) satisfies 

dtv = Cjv . 

We have u{t,x) = e*'-''uo{x) and v{t,y) — e'^'-'-'voiy) where uq — u(0, •) 
and vo = v{0,-) are related through vo{y) = ^/fpooiS^^^^y) uo{S^^/^y). In 
particular, it is easy to check that for alH > , 



UQl/joodx 



(det S) 



-1/4 



L2 



where 



(28) 



(^o(^'o))(y) = (2^)-^/^ f [ voiy)e-^^^'/' dy) 



-|yl74 
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is the L^-orthogonal projection of vq on the kernel Ce <> of £j. Thus, 
proving ([M)) is equivalent to proving 

e*^.^(/-77o)||^^^^^ < cf.{S)y'c.p (-^i) , (29) 
where here and in the following we use the standard operator norm 

ueL2(R«) ||w||l2 

for an arbitrary operator A. 

In the following, we will often work with Cj and Yt rather than with Cj 
and Xt- 



3 The linear algebra problem 

The stochastic differential equation ([S]) for the linear case (quadratic poten- 
tial) that we consider becomes 

dXt^-{I + J) SXt dt + V2 dWt , (30) 

and is associated with the drift matrix 

Bj:={I + J)S. (31) 

With the change of variables given in Section[5](Yt — S^^^Xt), the stochastic 
differential equation (150]) becomes 

dYt = ~{S + J)Yt dt + \/2 5^/^ dWt. 

The drift matrix is now 

Bj = S^f^BjS'^'^ ^{S + J), (32) 

where, we recaU, J = S^^^JS^^^ e AnW ■ We first collect basic spectral 
properties of Bj (or equivalently of Bj) when J e ^Ar(M)) and then show 
how this spectrum can be constructively optimized. 



3.1 Spectrum of Bj for a general J E A 



Proposition 1 For J e AnO^) (or equivalently J = S'^/'^JS-^/'^ e AnO^)) 
and S e (M) , the matrix Bj ^ S + J has the following properties: 

(i) a{Bj) C |z e C, Re(z) > o| . 

(ii) Tr(Bj) = Tr(5) . 
(iii) minRe[cr(B,7)] < ^p- . 
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Notice that the properties stated above on Bj also hold on Bj since (7{Bj) = 
a{Bj) and l^{Bj) = Tr{Bj) . 

Proof Let A G C be an eigenvalue of Bj with corresponding (non-zero) eigen- 
vector xx e C-'^: 

{S + J)x\ = Bjx\ = Xxx . 
Since S* is a real matrix, the complex scalar product with xx gives 

X\xxf = \S'/^xxf + ixx,Jxx)c. 

Here and in the following, the complex scalar product is taken to be right- 
linear and left-antilinear: for any X and Y in , 

Using the fact that J e ^jv(lR) , we get: 

This ends the proof of {i). The proof of (ii) follows immediately from the 
fact that the trace of the antisymmetric matrix J is . 
To prove (iii), let 

= {Ai,A2,...,A,} 

denote the spectrum of Bj , and let mk denote the algebraic multiplicity of 
Afe . Part (ii) says 



J2 mkXk = Tt{S) e 



and consequently: 



^mfcRe(Afc) = Tr(5). 
fe=i 

Now, using the fact that X^^^i mk = N , we conclude 

MS) 



min Re [(t(Bj)] = min {Re(Afc), fc e {l,...,r}} < 



N ■ 

□ 
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3.2 Optimization of cr{Bj) 

Our goal now is to maximize minRe[(T(-B,7)] over J G ^Ar(]R) , or equiva- 
lently, to maximize minRe[cr(i3j)] over J — S^^^JS^^^ € ^Ar(M) . Indeed, 
this is the quantity which will determine the exponential rate of convergence 
to equilibrium of the non-reversible dynamics (jl7p as it will become clear 
below. 

From Proposition[TJiii), the maximum is obviously achieved if there exists 
a matrix J S ^jv(]R) such that: 

yXe<j{Bj), Re(A) = ^. (33) 

In the following proposition we obtain a characterization of the antisymmetric 
matrices J (related to J through J — S^^^JS^^^) for which ((33|) is satisfied 
and Bj is diagonalizable (see ([M]) below). This characterization requires to 
introduce a companion real symmetric positive definite matrix Q e S^'^ (M) . 
The case of non-diagonalizable Bj is then discussed, using an asymptotic 
argument. We finally show how this characterization can be used to develop 
an algorithm for constructing a matrix J e An {^) such that ((55)) is satisfied. 

Proposition 2 Assume that J G .4Ar(M) and that S G 5^°(M) . Then the 
following conditions are equivalent: 

(i) The matrix Bj — S + J is diagonalizable (in <C) and the spectrum of Bj 
satisfies 

aiBj)c^^ + iR. (34) 

Tr (S^ 

(a) Bj — — / is similar to an anti-adjoint matrix. 

J" 

(Hi) There exists a hermitian positive definite matrix Q — Q such that 

JQ~QJ ^~QS- SQ+'^^^^Q. (35) 

(iv) There exists a real symmetric positive definite matrix Q = such that 
JQ-QJ ^-QS~ SQ+'^^^^Q. (36) 

Proof First we prove the equivalence between (i) and (m). Equation (j34p is 
equivalent to the statement that there exists a matrix P G GL„(C) (where 
GL„(C) denotes the set of complex valued invertible matrices) such that 



Bj —I ) P = diag(iti, . . . Hn) 



for some in M , which is equivalent to statement (m), since any anti-adjoint 
matrix can be diagonalized in C . 
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To prove that (ii) implies (iii), we write statement (ii) as: there exists a 
matrix P e G-L„(C) such that 

{^^y - ^/ = -P-'BjP + ^/ . (37) 
Since Bj = S + J e MnI^) and J G ^Ar(M) , we obtain 

p-^jp-p^j {p-'f ^ -p^s (p-'Y - P-^SP+^^I. 

— T 

We multiply this equation left and right by P and P respectively, to obtain 

JPP -PP J = -PP S-SPP + ^ ' PP . (38) 

X 

Statement {iii) follows now by taking Q = PP . Conversely, {Hi) {ii) 

follows from the writing Q = PP , with P e GL„(C) (take P — \/Q) for any 
hermitian positive definite matrix Q . Then, one obtains {ii) by going back 
from (l38)) to (j37p . Finally, {Hi) implies {iv) by taking the real part of ((35)) 
and using the fact that J and S are real matrices. The converse {iv) {Hi) 
is obvious. This ends the proof. □ 

Remark 5 Notice that if J is such that (IMl) is satisfied, so is — J (and thus 
J^). Indeed, if {J,Q) satisfies then (— J, Q^^) also satisfies ([55)) . □ 
Let us give another equivalent formulation of p6p. 



Lemma 2 With the notation of Proposition let us consider matrices J G 
•>°{R) and Q e S>° 



An{^), S e S^°{R) and Q g S^°{R). Let us denote {Xk}k=i the posit 



real eigenvalues of Q (counted with multiplicity), and {V'fcjfeLi the associ- 
ated eigenvectors, which form an orthonormal basis of M.^ . Equation ()36p is 
equivalent to the two conditions: for all k in {1, ... , N}, 

(^.,5V,)« = ^ (39) 
and, for all j ^ k in {1, . . . , N}, 

{Xj - Afc)(V'j-, JVfe)M = {Xk + Aj)(V'j, S'V'fc)R • (40) 

Proof Since {V'fejfcLi form an orthonormal basis of M.-'^ , Equation (|5S)) is 
equivalent to this same equation tested against ipj on the left, and ipk on the 
right. This yields: 

2Tr(5') 

XktpJJipk - Xj^pjjtpk = -XjtpJSipk - Xk^pJStpk H — SjkXk , 

where djk is the Kronecker symbol. When j = k , we obtain ()39p by using 
the antisymmetry of J , together with the fact that all eigenvalues of Q are 
non-zero. When j k , we obtain (|40)) . □ 
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Notice that when the eigenvalues of Q are aU with muhiphcity one, J is 
completely determined by for all j 7^ /c in {1. ... , N} , 

{tpj , JV'fe)R = + (v^j-,gVfc)M. (41) 

Afe — Aj 

Indeed, by the antisymmetry of J , the remaining entries are zero: 

(V'j, JV'j)k = for aU J e {l,...,iV}. 
This motivates the following definition. 

Definition 1 Wc will denote by VoptiS) the set of pairs {J,Q), where Q 
is a real symmetric positive definite matrix with N eigenvalues of multiplic- 
ity one and associated eigenvectors satisfying ([M)) . and J is the associated 
antisymmetric matrix defined by (j4ip . 

Notice that for any (J, Q) G VoptiS) , J is completely defined (by dH])) as 
soon as Q is chosen, so that the set VoptiS) can be indexed by the set of 
matrices Q S 5^°(M) with N eigenvalues of multiplicity one, and with eigen- 
vectors ipk satisfying ((39)) . As it will become clear below, the matrix Q of a 
pair (J, Q) £ VoptiS) appears in the quantitative estimates of Theorem[T]and 
Theorem [5] through the constants C^'' and C^"* . The construction of the pair 
(J, Q) is also better understood by splitting the two steps: (1) construction 
of Q and (2) when Q is fixed, construction of J . 

Remark 6 We would like to stress that the set VoptiS) does not provide all 
the matrices J € v47v(K) such that cr(i?j) C '^^jp-+iM. . Indeed, first, we have 
assumed that Bj is diagonalizable and, second, in this case we have assumed 
moreover that Q has N eigenvalues of multiplicity one. 

Actually the spectrum of Bj depends continuously on J. Hence any limit 
J = liuin^aoJn in ^Ar(M) with iJn,Qn) G VoptiS) wiU lead to aiBj) C 
'^^jf^ + iM. . A particular case is interesting: Fix the real orthonormal basis 
{i'j}jLi and consider Qa with the eigenvalues (a, . . . , a^) with a > . The 
unique associated antisymmetric matrix is given by iipj , Jaipj)m = and 

a'' + a-' 

i%pj , Ja1pk)R = z -ii^j , S''0fc)R • 

a'^ — 

Taking the limit as a — > -f 00 or a — > 0+ leads to 

ilp] , iooV'fc)R = -Sign(fc - j)(V'j , S'lpk)R, Jo+ = -Joo ■ 

Actually, for such a choice Jopt — Joo or Jopt = Jo+ j the matrix S -\- Jopt is 

triangular in the basis iiJ-'j)i<j<N and aiBj^^^) — | ^'^j^^ | • In general (see 

for example Subsection l3.3l) . the matrix Bj^^^ may not be diagonalizable over 
C and may have Jordan blocks. □ 
We end this section by providing a practical way to construct a couple 
(J, Q) satisfying (|36)) (or equivalently (J, Q) G VoptiS)), for a given S € 



17 



5^°(M). The strate gy is simple. We first build an orthonormal basis {'ipk}k=i 
of such that ([55]) is satisfied, then we choose the eigenvalues {Xk}^^i 
distinct and positive, and define J by (j41l) . The only non-trivial task is thus 
to build the orthonormal basis {ilJk}k=i ■ 

Proposition 3 For every S £ ^^"(M), there exists an orthonormal basis 
{i'k}k=i ofM.^ such that p9p is satisfied. 

Proof We proceed by induction on iV, using some Gram-Schmidt orthonor- 
malization process. The result is obvious for = 1 . For a positive integer 
N , let us assume it is true for — 1 and let us consider S G (R) . Let 
us set T = . The matrix T is in Sff°{R) with Tr(T) = 1 . Consequently 

(V'^,TV'.)E > 0, z = l,...iV and = 1 for any orthonormal 

basis {V'ilili of R^ . Assume that not all {ipi, Tipi)]Si are equal to . Then 
there exist iq, ii S {1,..., N} such that 

(V'zo > TiP.^ )k < ^ , ii^ti , T^p^, )r > ^ ■ 

Set ipt — cos{t)ipig + sin(t)-0.ij and consider the function f{t) ~ {^t,Tipt)M. ■ 
This function is continuous with /(O) < and /(7r/2) > Conse- 
quently, there exists a € (0, tt/2) such that 

iA,,T^tjM = ^. (42) 

Let now 77 = I — ipt, ("0*,)"^ denote the orthogonal projection to Spanl^ipt,)^ 
and define 

N - 1 

This operator is symmetric positive definite on Span(?/'4^)^ with 

Tr(Ti) = (Tr(T) - {^H„Ti;tM = 1 • 

It can thus be associated with a symmetric positive definite matrix in A4 Ar_i (R 
By the induction hypothesis there exists an orthonormal basis (V'2, ■ • • , '4'n) 
of Span(?/'4^)^ such that 

Let us consider the orthonormal basis of R^ : 
We obtain 



\ i> 2. 
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and, for i>2. 



This ends the induction argument. □ 

Remark 7 Finding t^, such that is satisfied yields a simple algebraic 
problem in two dimensions. Let (iq, ii) be the two indices introduced in the 
proof. The matrix , TV'jOr) .^^.^^^^^^^ j € X2(M) is 



ao P 



with ao < , «! > , /? € 



Then, G (0,7r/2) is given by 



(ai-^) (ao-^) 



and the vector ijjt, by 

i/'t. = , ^ =(V'»o + tani^Vii) ■ 
\/l + tan"' 

□ 

The above proof and Remark [7] yield a practical algorithm, in the spirit of 
the Gram-Schmidt procedure, to build an orthonormal basis satisfying (I39p . 
see Figure [T] This algorithm is used for the numerical experiments of Sec- 
tion [51 Notice that in the third step of the algorithm, only the vector ipn+i is 
concerned by the Gram-Schmidt procedure. The chosen vector ij^t, belongs 
to 'Sjpn ® ^'4'n+i and all the normalized vectors {tpn+2, ■ ■ ■ , "^Pn) are already 
orthogonal to this plan. 

A simple corollary of Proposition [3] is the following: 

Proposition 4 For every S G iS^'^(M) , it is possible to build a matrix J G 
y^Ar(M) such that 

= minRe[CT(i3j)] > minRe[cT(S')] 

where Bj = S + J . Moreover, this holds with a strict inequality as soon as 
S admits two different eigenvalues. 

In conclusion, the exponential rate of convergence may be improved by using 
a non-reversible perturbation, if and only if S is not proportional to the 
identity. We also refer to [23l Theorem 3.3] for another characterization of 
the strict inequality case. 
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Algorithm for constructing the optimal nonreversible perturbation 

Start from an arbitrary orthonormal basis {tpi, . . . ji/^jv). 
for n = 1 : iV - 1 do 

1. Make a permutation of (-(/)„, . . . , -i/jjv) so that 



and 



(^„,S'^„)m= max i^k,S^k)»>TriS)/N 

k — 'n,...,N 



{il)„+i,Stpn+i)R = min (-i/jfc, Si/)fc)M < Tt{S)/N . 



2. Compute t, such that tpt, = cos{t,)tp„ + sin(t*)7/;„+i satisfies (V't, , 5'i/'t. )m ~ 
Tr{S)/N (see Remark [7] above). 

3. Use a Gram-Schmidt procedure to change the set of vectors {tpt, , V'n+ii ■ • • > '4'n) 
to an orthonormal basis {tpt^ , ipn+i, • • • , i^n) ■ 

end 



Fig. 1 Algorithm for constructing the optimal nonreversible perturbation 



3.3 Explicit computations in the two dimensional case 

In the two dimensional case {N — 2), all the matrices J such that (t{Bj) C 
Ti{S)/N + iM. can be characterized. Accordingly, explicit accurate estimate 
of the exponential decay are available for the two-dimensional ordinary dif- 
ferential equation: 

^ = -(/ + J)Sxt with xo given in . (43) 

After making the connection with our general construction of the optimal 
matrices J (see Definition [ij , we investigate, for a given matrix S £ 5^°(M) , 
the minimization, with respect to J , of the prefactor in the exponential decay 
law. We would like in particular to discuss the optimization of the constant 
factor in front of exp(-Tr(S')t/2) . 
Without loss of generality, we may assume that 



S = 



1 
A 



and J = 



a 

-a 



where A > is fixed. The eigenvalues oiBj belong to Tr(S')/2 = (l-t-A)/2-f-iM 
if and only if 

and then, the eigenvalues of Bj = {I + J)S are 



A + 1 ± iJAXa^ - (1 - A)2 

'^± = 2 ■ 

When the inequality (|44p is strict, the associated eigenvectors are 



1 \ /i±-l A-l±zV4Aa2-(l-A)2 
u± ^ { with a± = ; — = — — 
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P- 



The matrix Bj equals 

B.,=P 
with P = 

The case a — zfc gives the matrix 

Bj = 



1 1 



and 



a_ -1 
— q;_I_ 1 




which has a Jordan block when A 1 . This ends the characterization of all 
the possible optimal J's in terms of the exponential rate. 

Let us compare with the general construction ofthepair((3,J = JS^), 
see Definition [T] The matrix Q is diagonal in an orthonormal basis (V'ijV'2) 
which satisfies the relation ([55]) . This yields 



2|2 



1.1 12 




Then, 



Up to trivial symmetries one can fix t/;! = [ ^ | and '4'2 
from (Uni), the eigenvalues Ai, A2 of Q must satisfy 

(A2 - Xi){-2aVX) = (A2 + Ai)(l - A) 

and the limiting cases a = are achieved only after taking the limit 

^ ^ +00 or 

and consider the two-dimensional Cauchy prob- 



Assume now > 



-00 . 

4A 



lem Its solution equals 

exp(-^+t) 



Xt= P 









exp(— 



which leads to 



l^t||<ll^llll^-'l 



exp 



1 + A 



lla^oll, 



when II • II denotes either the Euclidean norm on vectors or the associ- 
ated matrix norm, ||yl|| = ^J\nsxia\A*A)) . This yields the exponential 
convergence with rate Tr(S')/2 = (1 + A)/2, as soon as a satisfies > 



4A 



, while the degenerate case a 



2 _ (l-A)^ 
4A 



would give an upper bound 



C{l~\-t)e 2 * . A more convenient matrix norm is the Frobenius norm given 
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by \\A\\p = J2i.j=i\^i]? J2ae<7{A'A)^ '^i^h the equivalence in dimen- 
sion 2 , ^11 All J. < II All < II All F. By recalling a+ = — , we get 

||x*||<||P||i=^||P-^||Fexp (-^i) ll^oll 

„(l + |a+P) / 1 + A \ „ „ 
<2^ ^^exp :^t] xo 

Now, it is clear that the infimum of ||P||i?||P^^||i? is obtained asymptotically 
as |a| — > oo and equals . It corresponds to an antisymmetric matrix J 
with infinite norm. 

To end this section, we would like to discuss the situation when the orig- 
inal dynamics (when J = 0) has two separated time scales, namely A is 
very large or very small. In the case A 1 , we observe that the optimal 
II-P||f||^'"^||f (and thus the optimal ||P||||P"D scales like , and that 

this scaling in A is already achieved by taking ~ ^"'"^^^ (twice the mini- 
mum value in (04])), since in this case, ||P||f||^' "'^||_f = V2- — ^ . In terms 

vA 

of rate of convergence to equilibrium, it means that, to get ||a;t|| of the order 
of ||a;o||/2, say, it takes a time of order ln(l/A). This should be compared to 
the original dynamics (for a = 0), for which this time is of order 1/A. Of 
course, a similar reasoning holds for A 3> 1 . Using an antisymmetric per- 
turbation of the original dynamics, we are able to dramatically accelerate 
convergence to equilibrium. 



4 Convergence to equilibrium for Gaussian laws and applications 

In this section, we use the results of the previous section in order to under- 
stand the longtime behavior of the mean and the covariance of Xt solution 
to dSOl): 

dXt = -{! + J)SXt dt + V2dWt. 

In particular, if Xq is a Gaussian random variable (including the case where 
Xq is deterministic), then Xt remains a Gaussian random variable for all 
times, and understanding the longtime behavior of the mean E(Xt) and the 
covariance matrix Var(Xt) = E(Xt (g) Xt) — lE(Xt) (gi lE(Xt) is equivalent to 
understanding the longtime behavior of the density of the process Xt, which 
is exactly Corollary [T] in a very specific case. Here and in the following, (g) 
denote the tensor product: for two vectors x and y in M^, x ® y — xy^ is a 
N X N matrix with (i, j)-component Xiyj . 
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4.1 The mean 

Let us denote Xt = E(Xt) , which is the solution to the ordinary differential 
equation 

^=-{I+J)Sxt, xo = x. (45) 

The longtime behavior of xt amounts to getting appropriate bounds on the 
semigroup e~i^+J)St equivalently on e~('^+'^)* . 

When J = , namely for the ordinary differential equation 

dXt r, 

— =-Sxt, xo = x, 
at 

we immediately deduce from the spectral representation of the positive sym- 
metric matrix S that 

II a;* II < e-''*||a;o|| , (46) 

where 

p := min{cr(5)} . 
The above bound implies that 

||e-^*|| <e-''*, 

where ||M|| = sup^gKJv^^^o ■ Notice that p < ^p- . 

We now derive a similar estimate for the semigroup generated by the 
perturbed matrix Bj = S+J (or equivalently Bj — {I+J)S), when (J, Q) € 
Vopt , and show that a better exponential rate of convergence is obtained. As 
explained in the introduction, the price to pay for the improvement in the 
rate of convergence is the worsening of the constant (which is simply 1 in the 
reversible case) in front of the exponential. 

Proposition 5 For {J,Q) G Vopt o,nd J = 5~^/^J5~^/^, the estimates 

-{S+J)t 



e 



^-{I+J)St 



<K(Q)V2exp(-^i) , (47) 



<«(Q-i5)i/2exp(-^t 



(48) 



hold for every t>0. 

Proof Consider the ordinary differential equation 



^ = -{S + J)yt, yo = y. (49) 
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We introduce the scalar product (•, ■)q-i := (-,(3 ^■)r on with the cor- 
responding norm || • ||q-i . We calculate: 

j^\\yt\\l-i ^ ~2{yuQ~\S + J)yt)m 

= - {yt, {Q-^S + Q-^J + SQ-^ - JQ-^) yt)^ 
2Tr(5) ^_i 

yt, — — Q yt 

2Tr(5)„ ,,2 
= ^llytllQ-i- 

In the above, we have used the identity 

2Tr(S') 



N 



which follows from after multiplication on the left and on the right by 
Q^^. From the above we conclude that 



We now use the definition of the norm || • ||q-i to deduce that 

htW < WQ'/'WWQ-'^'ytW 
= \\Q'^'\\\\yt\\Q-^ 
<e-^*||Q^/^||||y|lg-. 

<e-^*||QV2||||Q-i/2||||y|| 
= e-^*^(Q)i/2||y||. 

For the second estimate, we set xt ~ S^^l'^yt and obtain 

|k,||<||5-V2Ql/2||||g-l/2j^^|, 

<e-"-^'K{Q-^SYl^\\x\\. 



□ 

Proposition [S] shows that, for a well chosen matrix J, the mean xt — ^(Xt) 
converges to zero exponentially fast with a rate "^^^ ■ Equation ([2^ in 
Theorem [T] is a simple corollary of and the inequality k{Q^^ SY^^ < 
K{Qy^^K{S)^^^ , so that C^^^ in (gH) can be chosen as 

C^r^^^KiQy/^ (50) 
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Remark 8 Let us make a remark concerning the constant C^' in (|22[) . using 
the upper bound (|50p . It is possible to have Cn independent of N , while 
keeping the norm of the perturbation J under control. More precisely, for a 
given orthonormal basis {ipk) satisfying (j39p . let us consider the eigenvalues 

Xk ~ N + k . On the one hand, C^'* remains small since k{Q) — 2 . On the 
other hand, using (|^T|) . we have 

j<k ^ ■' -' 

<2(47V)2^(V'„^Vfe)R 
j<k 

< 16N^\\S\\1. 

Thus, the norm of J (compared to the one of S) remains linear in . □ 
4.2 The covariance 

Let us again consider Xt solution to ((30|) . and let us introduce the covariance 

St = ® Xt) - EiXt) ® E{Xt) , 
which satisfies the ordinary differential equation: 

-^^^iI + J)SUt-StSiI-J) + 2I. (51) 

The equilibrium variance is Soo — S^^ . 

Proposition 6 For {J,Q) £ Vopt and J — S^^^^JS^^^^ , the estimate 

||^t-5'-^|| < K(Q-i5)exp (^-2-^^i^ \\Eo - S'^\\ (52) 

holds for all t > , when the matricial norm is induced by the Euclidean 
norm on . 

Proof The solution to ^ (see e.g. [1H1I3S]), is 

St = S-^ +e-'^''iSo - S-^)e-'^^. (53) 

The result then follows from the estimate on ||e^*^'|| in Proposition [5] above 
and ||e-*^^|i = ||(e-*^'')^|| = \\e-'^-'\\ . □ 



25 



4.3 Gaussian densities 

As a corollary of Proposition [S] and Proposition [HI we get the following con- 
vergence to the gaussian density (see (jl6|) ) 



= — (2^)JV/2 — ''^P 1 ^ ) ' ^^^^ ^oo ^S. 



Proposition 7 Assume that Xt solves (1301) while Xq is a Gaussian random 
variable, so that Xt is a Gaussian random variable for all time t > , with 
the density tp^ . Assume moreover that J — S^^^^JS^^^^ , and that {J,Q) 
are chosen in Vopt ■ Then, the inequality 



< N2^e 



'exp(2e-2^(*-*«)||xo||2 



holds for all times t larger than 



in = 



N 
2Tr5 



In [A{l + \\S\\)^{Q-^S){l + \\SS4)\ 



(54) 



This result is related to the result stated in Corollary [1] that will be proven 
in Section [S] Corollary [T] provides a better and uniform in time quantitative 
information (which has also a better behavior with respect to the dimension 
N according to ((68)) and Remark ITOl) . On the contrary, it requires more reg- 
ularity than Proposition [7] which does not assume V'o ^ L'^(^^ t'^^ dx) . Of 
course, with initial data outside i^(M^, -0^^ dx) , the convergence estimate 
makes sense only for sufficiently large times (hence the introduction of the 
positive time t^ in Proposition [7]). 

Proof The Gaussian random vector Xt has the mean xt , which solves (pSj), 
and the covariance St , solution to (ISTI) , so that 



det(Z't 



N-l/2 



(27r)^/2 



■ exp 



(a; - xt)'^Ut ^{x - Xt) 



When t >to , Proposition|6]gives HiTt— i7oo|l < 4 H^^'ooH and thus, \\Z!co'^ StSoo'^ - 
/|| < i , which yields 1^00 < St < and < S^' < . In 

particular < allows to compute 



l + ||V'/-^/'oo|l'..^-U =1+ / (V'/-V'oo)'fc' = 



-N/2 



det(i:t 



det(Z'oo)^^/^ Jrn 



exp 



{x - Xt)'^St ^(^ ^ ^t) + 



(A 



J\2 



T V— 1 ^ 

X ^ X 
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We then use the relation, for A and B in 



{x — Xt)'^A(x — xt) — x'^ Bx 

= {x-{I- A-^B)-^xtf{A - B){x A-^B)-^xt) 
+ xj [A- A{A- B)-^A] Xt, 

with A = and B — in order to get 



\ -1/2 



X exp [ 



Xt 



detiS^^St)^ det(2/ - S^^Et)^ 

X exp [2(2/ - U^'St)'' - l] ^-^Xt) 



After setting Rt — I — S^St , we deduce 

IIV'/-V'oo|li.(^-l) 

1 



det(/-i?2)i 
1 



det(/-i?2)2 



r X [exp {xj[2{l + Rt)-^ - I])S^^xt) - l] . (55) 



Let us start with the determinant det(/ — i?j ) . the condition t > to and 
Proposition [6] give 



\\Rt\\ = \\I-E^'St\\< 
With II < , we know 



— and ||i??||<- 



-4^(t-to) 
16 



16 

I hidet(/ -Rf)\< |Tr(ln(/ - Rf))\ < N\\ hi(/ - i?2)|| < -7Vin(l - ||i?2||) 
We deduce 

1 



det(/-i?2)2 



-<(l-||i?,1|)-^< 



< 2~ . 



Concerning the exponential term in (|55p . Proposition [5] implies that the ab- 
solute value \xf [2(7 + _Rt)-i — /] Z'j^^xtl is smaller than 



(l + 2||(/ + i?)-i||)||i:r^|| x«(Q-i5)exp 



-2^t)||.o|P 
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The inequality ||(1 + R)-^\\ < (1 - < | and the condition t > to 

imply WSr'W < lll^oo'll = lll^ll and 

\xj [2(1 + I] S^'x,\ < -ii-e-2^(*-*o)||:,o|i2 



9x4 



We have proved 



llV't'-^'ooll 



-u < 



, Tr(S) 



(t-to) 



2- 



16 / 

11 nTr(S) 



exp 1 -e 



By using (1 - a;)-^/^ < 1 + N2^/'^x when x € (0, 1/2) for the first term, and 
— 1 < ye^ when y > for the second term we finally obtain 



« -4^(t-to) 



11 „ N n Tr(S) , > 

-_-2"2"e~ "^1 



11 

a^oirexp [ —e 



-(*"*"^iixoir 



which yields the result. 



□ 



4.4 General initial densities 



As a corollary of Proposition [71 a convergence result for a general initial 
probability law can be proven by using an argument based on the conditioning 
by the initial data. 

Proposition 8 Let ■0/ satisfy the Fokker-Planck equation (j25p . with an ini- 
tial probability law with density and such that f-^f^ e""^" '0o '^•^ ^ 
/or some positive a . Assume moreover that J = S^^^^JS^^^^ , that {J,Q) 
are chosen in Vopt o,nd that to is given by (j54p . Then the inequality 

II^Z-^o^ll^ <7V2^+V^^(*-*") / (56) 

holds for all t > ta — to + 2Tr{s) I )l ■ 

Proof In all the proof, J = S~^^^JS~^^^ is fixed, with (J, Q) chosen in Vopt ■ 
For X e and t > 0, let us denote 4>f the density of the Gaussian process 
X'^ solution to: 

dX'f = -(/ + J)SX^ dt + V2dWt with X^ ^x. 

Proposition [7] can be applied with = 0f and i^o = , so that the time 
to = 2^ In [4(1 + ||S'||)k(Q~^S')] is fixed. With the decomposition 

i^tiy)^ f rt{y)t^{x)dx, 
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coming from (j)Q — Sx , we can write: 



IIV'Z-V'oo II' 2(^-1) 



1 



< 



f my)? 
R« Jrn ipoo {y) 



dy-l 

<t>t{y)to{x)dx^ dy-l 
dy ■00 (^) dx — 1 



dy — l\ -0Q (x) dx . 



With Proposition [71 we deduce 



ll0/-0oo|li.(^-i)<iV2"e-^^ 



1 



|x|| exp 2e 



N 



1 + - / e""^" da; 



for t > to = + 2T>(g) |ln(f )|- To get the second hne, we used (for t > 
ta) e^2^(t-to) < I and (for u > 0) uet" < ^e"" < -ie"" . Writing 
e~ = e~ "iv^*- °''e"l "tI and discussing the two cases a > 4 

and a < 4 yield the result. □ 
The aim of the analysis using Wick calculus in Section [5] is to obtain more 
accurate and uniform in time estimates. 



5 Convergence to equilibrium for initial data in L^{M.^ ,xpaodx) 

We shall study the spectral properties, and the norm estimates of the corre- 
sponding semigroup, for the generator Cj defined by (I27|) (with y replaced 
by x as a dummy variable): 

£j = V'^SV - ^x^Sx + ^Tr(S') + i(a;'^Jv - V^Ja;) . 

The operator Cj acts in L'^{M.^ , dx;C) and is unitarily equivalent (when 
J = S^^^JS^/^, and after a change of variables, see Section to 

Cj = -{BjxfV + A with Bj = {I + J)S, J^An, 

acting on L^(M^, ^Aooda;; C) . Since for J ^ , the operator Cj (or Cj) is 
not self-adjoint, it is known (see [40 . 8".' T4llT9l[T7l[T8] ) that the information 
about the spectrum is a first step in estimating the exponential decay of the 
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semigroup, but that it has to be completed by estimates on the norm of the 
resolvent. This will be carried out by using a weighted L^-norm associated 
with the constructions of the matrices Q and J introduced in Section |31 



5.1 Additional notation and basic properties of the semigroup e^^-' 
Let us introduce some additional notation. 

— We choose the right-linear and left-antilinear convention for L^-scalar 
products (or S — iS'-duality products) : 

(/,ff)L2 = / fix)g{x) dx. 

— For a multi-index n — (rii, . . . , njv) G , we will denote n\ = Y[f=i "-i^ 
\n\ — X^jLi "-i ^^"^ when Xi, . . . ,Xn belong to a commutative algebra 

— The space of rapidly decaying complex valued C°° functions is 

sup \x"d^f{x)\ < C^p] 

and its dual is denoted S'{M^) . 

— The Weyl-quantization q^{x,D.^) of a symbol q{x,^) G ^'(M^^) is an 
operator defined by its Schwartz-kernel 

For example, for q{x,£,) = f{x), {x, D^) is the multiplication by f{x), 
for q{x,£_) — /(C), q^{x,Dx) is the convolution operator f{—i\7), and 
for qix, = q'^ix, D^) is ^(x^V + V^x) . 

— The Wick-quantization of a polynomial symbols of the variables (z,z), 
where z £ is an operator defined by replacing Zj with the so-called 
annihilation operator aj = dx^ + ^ and 'Zj with the so-called creation 
operator a* = —dxj + ^ . Wick's rule implies that for monomials involving 
both z and z, the annihilation operators are gathered on the right-hand 
side and the creation operators on the left-hand side: For given multi- 
indices a, (3 € , the monomial z"z^ becomes (a*)"a^ . The properties 
of the Wick calculus that we need here are reviewed in Appendix [X] We 
shall also use the vectorial notation 











a — 




, a = 






\aNj 







with their transpose and a*''^ . 
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— The orthogonal projection from L^(M^, dx; C) onto Ce « will be de- 
noted by Uq . 

Let us now recall a few basic properties of the semigroup e*'^-'. The Weyl 
symbol of 

-jCj + = -V^5V + Ix'^Sx - lix^JV - V^Jx) 

2 4 2 

is (using the fact that J is antisymmetric) 

= es^ + ^ - ^ ix^j^ - ejx) = es^ + ^ - ^x^J^ , (57) 

which is a complex quadratic form on . Besides, the operator — £j is the 
Wick quantization of a quadratic polynomial since 

-jCj = -V^SV + ^x^Sx - ^Tr(5) - ^(a;^^V - V^Jx) 

= a*'^{S-J)a. (58) 



Proposition 9 The differential operator —Cj is continuous from S{M.^) 
into itself and from iS'(]R^) into itself. Its formal adjoint is —C-j ■ With the 
domain D{-Cj) = [u & L'^{W^ , dx; C), -Cju e L^R^, dx; C)} , the opera- 
tor —jCj is a maximal accretive and sectorial operator in L^{W'^ ,dx;C) . Its 
resolvent is com,pact and its kernel equals Ce * . The associated semigroup 
{e*^-')t>o has the following properties: 

— For any u e <S(M") (resp. any u e S'iW^)), the map [0, +00) ^ t ^ 
e^'^^u is a S{R'^) -valued (resp. S'{R^) -valued) C°° function. 

— For any t > 0, the operator e^'~-^ sends continuously S' {M.^) into S(M.'^) . 

— In the orthogonal decomposition L'^{M.^ ,dx;C) = ^^en'^k the finite 
dimensional vector spaces spanned by Hermite functions with degree k: 

i(a*)V-^ --^^ 



Vk = Span ^{a* , n G N"^ ,\n\ = k} , 

the semigroup has a block diagonal decomposition 

± 

„tCj ^ ^ I 



e 



Proof As a differential operator with a polynomial Weyl symbol, —Cj is 
continuous from 5(]R^) (resp. 5'(]R^)) into itself. Its formal adjoints has the 
Weyl symbols qj{x,£^) = q-j{x,S,) and equals —C-j. For k &N, set 

= {u e L'^{R^, dx; C),x"D^u e L'^{R'^, dx; C), 

for all a, /3 e s.t. \a\ + |/3| < fcj 
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and let Ti he its dual space. They satisfy 

n H'' = 5(M^) and U U'' = 5'(M^) . 

fcez feez 

Since 5* is a real symmetric positive definite matrix, the inequality 

\qj{^,o\>es^+^>cs{\e + M') 

implies that the operator —Cj is globally elliptic (see fTBl[551[551l5B] ) . There- 
fore, it is a bijection from H'^ onto 'H''~^ for any fc e Z. This provides the 
compactness of the resolvent and the maximality property. The sectoriality 
(see [371 Chapter VIII]) comes from 



(w , —Cju)]^2 — {u , a*'"^ Sau) — {u , a*'^Jau) 
with \{u , a*'"^ Jau) ]^2 \ < — ^ — 777\(''^' a*''^ Sau) 1^2 . 



L2 



min <j{S) 

This yields (using the fact that Re((u, -~Cju)i^2) — {u , a*''^ Sau) [^2 and 
Im((u, -Cju)l2) = -(u, a*^'^Jau)L2) 

Vu e 5(M^), |arg (u , ~Cju)l2 \ < 6* , (59) 

with < tan(0) < „,ii['^^[g-) < +00. Here and in the following, arg(z) denotes 
the argument of a complex number z . 

Then the usual contour integration technique for sectorial operators (see for 
example [371 Theorem X.52] and its two corollaries) implies that {~Cj)^e*^-' 
is bounded for any fc e N and any t > Q . Combined with the global ellipticity 
of —Lj , this provides all our regularity results. 

The orthogonal decomposition L^(M^, dx; C) = ©^^pj 2?fe is actually the 
spectral decomposition for the harmonic oscillator Hamiltonian a*'^ a. From 
the Wick calculus (use either [ai,a*] = 5ij or the general formula recalled in 
Proposition 1131 3 in the appendix), we deduce 

[a*'^a , a*^^{S - J)a] = a*'^ [l AS ~ J)] a = 0. 

This implies that the spectral subspaces Pfc, fc g N, are indeed invariant by 
the semigroup e*'-' ' . □ 
Note that with the last property, the question of estimating the conver- 
gence to equilibrium stated in Theorem [5] is equivalent to estimating the 
decay of the semigroup 

e*^^(/-i7o) or e*^'^L^ 

I a; I 2 

where ilo is the orthogonal projection onto Ce ~ = 2?o (see also ((29)) '). 
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5.2 Spectrum of Cj 



The result of this section is a direct apphcation of the general results of [20l 
developed after [38',l^. See also 134,33 where these general results 
are used in order to compute the spectrum of the generator of a linear SDE 
with, possibly degenerate diffusion matrix. This result was first obtained 
in [30] using different techniques. 

Proposition 10 The spectrum of the operator ~Cj equals 




and its kernel is Ce -i . 

Proof The spectrum of the operator qY{x,Dx) — —Cj + -^^^ associated 
with the elliptic quadratic Weyl symbol qj{x,£_) defined by (|57)) equals, ac- 
cording to [20 -Th 1.2.2, 



J2 -iX{rx + 2kx),kxe 



A e cr{F) 
ImA > 



where F is the so-called Hamilton map associated with qj, and r\ is the 
algebraic multiplicity of A G o-{F) , i.e. the dimension of the characteristic 
space. The Hamilton map is the C-linear map F : — ^ associated 
with the matrix 

^ / 

-/ 



F = 



M, 



where 



-- / 



iJ S 



e M2n{C) 



is the matrix of the C-bilinear form associated with qj . The matrix F is 
similar to F defined by 



F 



Thus, the characteristic polynomial of F can be computed by 



r 1 


" 


F 


'V2 " 


1 


'iJ S' 





V2. 




^ 2 


_-S iJ_ 



det(F - A/) = dct(F - A/) = 2 

= 2"2Ar 



-2N 



iJ - 2X1 S 
-S iJ-2XI 



iJ - 2X1 S 
-S -J- i2XI i{J + S + i2XI) 

i{J-S + i2Xr) S 

i{J + S + i2XI) 

= 2-^" det(5 - J - i2XI) det(5 + J + i2XI) 

det(5 + J - i2XI) det(5 + J + i2XI) 



-2N 

2 
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where we used det(Af ) = det(M^) for M = S - J - i2XI in the last Une. 
Using the fact that Re{a{Bj)) > 0, we thus obtain that a{F) n {A, ImA > 0} 
equals ^(j{S + J) — ^a{Bj) . In particular one gets, 

ImA > A'^- \ J ^ 

and 

Tr(^) 



ImA > 



This concludes the proof. □ 
The Gearhart-Priiss theorem (see [TBHT^HU] ) provides the following corollary. 



Corollary 2 When the pair {J,Q) belongs to Vopt, the spectrum of —Cj is 
contained in 



{0} U <^ z e C, Rez > 



and 



where, we recall, 



lim In 



N j 

Tr(^) 



e*^''(/-i7o) 



N 



The above logarithmic convergence is weaker than an estimate ||e*'''-'|| < 

Tr(S) 

Ce «— with a good control of the constant C. Obtaining such a control is 
not an easy task for general semigroups with non self-adjoint generators (see 
[T51I191I17II14] ). This is the subject of the next section. 



5.3 Convergence to equilibrium for e*' 

Consider a pair (J, Q) e Vopt according to Definition [T] We recall that 
{J,Q) € Vopt satisfies ((M)) . With the matrix Q, we associate the operator 

CQ^a*'^Qa, (60) 

with which a natural functional space will be introduced in order to study the 
norm of e*'^-' . The operator Cq is the Wick-quantization of the polynomial 
z'^Qz. 

This operator Cq has the following properties: 

— It is continuous from 5(M^) into itself and from iS'(M^) into itself. 

— It is globally elliptic (see and it has a compact resolvent. 

— It is a non negative self-adjoint operator in L^(K^, dx; C) with the domain 
D{Cq) = {ue L^(R'^,dx;C),CQU e L^iR'^ ,dx;C)} . 
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— Its kernel is Ce ^ . 

— It is block diagonal in the decomposition L^(M^, dx; C) = 2?^-: 



feeN 



(61) 



One defines the two Hilbert spaces: 

— Hq = {w e L^(R^ ,dx;C), (u , Cqujl'i < +00} , naturally endowed with 
the scalar product 



L2 



- H^Q = W'qDV^ (where, we recall, = C 
product 



4 ) endowed with the scalar 



Proposition 11 Assume that the pair (J, Q) belongs to Vopt ■ Then the 
semigroup {e'^^-')t>Q is a contraction semigroup on T-Lq satisfying the fol- 
lowing estimate: 



Vi > 0. 



e*^'(/-i7o 



< e « 



(62) 



whe 



e'^'il^no] 



= sup p^i 



Proof The operator e*'''J is block diagonal (see Proposition IH]) in the decom- 
position Dq Dq ~ ©fegN -^fc which is an orthogonal decomposition in 
L'^iR^,dx;C) and also in owing to With e'^'e"-^ = e"^ , the 
semigroup property on is thus a consequence of the estimate (j62p in . 
Using the relation ((36|) together with the inequality ([76]) of Lemma |4] in the 

Appendix, we have: for all u G P = C[a;i, . . . , XN]e ~ H "Dq , 

2Tr(S') 

{u , {-^jCq - CqCj) > — {u , Cqu)l2 . 



Since the semigroup {e^ •')t>o is a strongly semigroup on 5(ffi^) and 
leaves V C 5(M^) invariant, we can compute for any u G V , 

d 



The proof is then completed using the density of D in Tig . 
We are now in position to state the main result of this section. 



□ 
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Proposition 12 Assume that the pair (J, Q) belongs to Vopt ■ Then the 
semigroup {e*^'')t>^ satisfies: 



Vi > 0. 



max(T(Q) 



minA,A'eCT(Q) , A#A' |A - A' 

Proof From the inequalities on real symmetric matricesmin (t((3) / < Q < 
max(T((5)/ and minCT(5)/ < S < maxa{S) I, we deduce with the help of 
Proposition [T3]-l the following inequalities on self-adjoint operators 



and 



mma{Q)a*'^ a < Cq < maxcr((3)a*'^a , 

na,xa{S) 

max (t((5) Cj + Cj 



™- /'ox *.T ^ ^ ^ I o\ *.T 

mma[b)a ' a < < ma,xa(b)a ' a 



mincr((5) Cj + Cj 
maxa{S) 2 



<Co< 



min (t(S') 2 



Here, we have used the fact that — is the Wick quantization of (z, S'z)c 

(see Proposition [T51 and Lemma S] below). Hence, using Proposition [TTl the 
following inequalities hold: for any u £ Wq and any < > io > 0, 



max(T(S') 
min(T(Q) 

2 



< niaxa(^) ^_,T4S)(^_^^^ 
~ min (t{Q) 

, , TTr(S), „Tr(S). t r C, 1 A" C*j t r 



L2 



L2 



Using the inequalities 



VweX>^, mincr(S')||u|||2 < (w, ^-^^ — -v) <\v\i^i\L,jv\^^ , 



L2 



with V — e'^'-'-'u and w = e*°^-'u, we deduce 



n Tr(S) ^ 
Tr(S) J e « 



t2mincr(S')2 



L2 



By taking to = ^ maxV(s) ' '^'^ obtain, for all ueHg 



\e'^'u\\i2 < K{Q)K{SYe-'—^e' sup 

t'>0 



£(L2) 



" r2 ■ 



(63) 
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The Lemma [3] below provides the bound 
with a e (0, 7r/4) defined by 



sup 
t'>o 



2 1 

< (64) 

TT^sin'^a 



, \ |Im {u,Lju)l^.\ \\J\\ 
tan — — 2a = sup < 



u<^D{Cj) |Re {u,Cju)l2\ mincr(S') 
The last inequality was proven in (ICTl) above. We thus obtain 



1 ^ 2cosa lUI 



sin a cos(2q;) min (t(S') 

In view of one can assume that a E (0,7r/8) (up to changing a by 

min(Q;, 7r/8)) so that 

^ < 2n/2 . (65) 



sin a min(T(5) 

When {J,Q) E Vopt (see Definition [1]), the relation (PT|) provides an ex- 
pression of the linear mapping associated with J in the orthonormal basis 
{4'k)i<k<N ■ In this basis, the Frobenius norm can be computed and we get 



\'<\\J\\% 



< 2 



2 

max(T((5) \ ||„||2 



lA-A" ' ""''^ 



i^)^max(.(^))= 



<2(- "'"''''"^^^ -] Nnmx{a{S))\ (66) 

,mmA,A'ecr(Q) , A#A' 



By gathering (l63|) - (|64)) - ((65)) ~ (l66t . we finally obtain the expected upper bound 
when t >to ■ 



||e*^-^u|li. < 2i"iV2K(Q) 



maxcr((3) 



minA,A'eff(Q).A#A' |A- A'| 
xniSye-'^'\\u\\l., 

for all u G and by density for all u G Vq . When t < to — 7^ , simply 
use 

|le*^^(/-7To)|U(L2) < 1 

< 2^N.iQ)i (- -^---(Q) )\^s)i X e-. 

VmmA,A'Ga(g),A5^A' |A- A'ly 

□ 
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Remark 9 A lower bound can be given for || J|| with 

\\jr>^\\j\\F 

min a{Q) \ 3 

- N Vnmx,,ve.(Q),A^A'|A-At' H 

min (t{Q) 



N VmaxA,A'6fT(Q),A5<tA' |A- A'l 



Tr(5'). 



Thus, we have 



J ^ min o-(g) min cr(5') 

maxAa'e£7(Q) , a^^A' |A - A'| ' 



(67) 
□ 



Lemma 3 Let (L,D{L)) be a maximal accretive and sectorial operator in a 
Hilbert space % with 

Vu e D{L) , |arg(u , Lujnl <S= — -2a with a > , 

where, we recall, arg(z) denotes the argument of a complex number z . Then, 
the associated semigroup satisfies 



Proof The case t = is obvious. 

For t > 0, e~^^ sends % into D{L) so that tLe~*^ belongs to . Consider 
first the case when ^ cr(i) . Our assumptions with a > , ensure that the 
operator tLe~*^ is given by the convergent contour integral 

tLe-*^ = — [ tze-^'^iz- L)-^ dz, 

where F is the union of the two half lines with arguments f — a and — ^ + a . 
For z = X ± G r with x > the resolvent (z — L)~^ satisfies (see 

for example [371 Chapter VIII.17]) \\{z - Ly^Wccu) < i • Moreover, \dz\ = 

1 + , \ dx = and |e~*^| = e^*^ . From these estimates, we deduce 

tan^ a sm a ' ' ^ 

" 27r ,/n sma a;sma Trsin a 



TT sin^ a 



When G cr(L) it suffices to replace L by e + L which satisfies the same 
assumptions as L with the same a with ^ a{e + L) . The identity 



t{e + i)e~*(^+^) - Le-'^ = tee-'*e-*^ + (e''* - l)tLe' 



tL 
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with t > fixed and e"*^ , tLe^*^ G CCH) implies lim^^o ||t(£ + i)e"*(''+^) - 
iie"*^||£(-H) = , which yields the result in the general case. □ 
In view of (|29p , Proposition [T^] yields the estimate (|24l) in Theorem [2] 
with a constant 

= 2^iV.(g) ( _ -^---iQ) \ ' . (68) 

To conclude, let us comment on the way of Cat behaves. 

Remark 10 In view of the upper bound (j68p . using the same construction as 
in Remark ISl we again notice that it is possible to have — 0{N'^) while 
keeping a reasonable perturbation J (with a Frobenius norm estimated by 
\\J\\f < 4A^||S'||_f). Contrary to the case of the ordinary differential equation 
discussed in Remark[5]our estimate does not provide a uniform in N constant. 

□ 



6 Numerical Experiments 

The algorithm for obtaining the optimal non-reversible is presented as a 
pseudo-code in Figure [T] 

In this section we present some numerical experiments, based on the al- 
gorithm presented in Figure [T] The numerical computations presented in this 
section are based on the following steps: 

1. Calculate the orthonormal basis {^k}k=i using the algorithm presented 
in Figure [T] 

2. Choose the eigenvalues of the matrix Q, {Xk}^^i, e.g. according to Re- 
mark [S] 

3. Calculate the optimal perturbation J using (|4T|) and the formula J = 

4. Calculate the optimally perturbed matrix Bj — {I + J)S. 

5. Calculate the matrix exponentials e^*^ and e"*^-' and their norms. 

In Figure [2] we present the results for a two dimensional problem, for 
which all results can be performed analytically, see Section 15751 We consider 
the case where the matrix B has a spectral gap, 

S' = diag(l,0.1). (69) 

In the figure we plot the norms of the matrix exponentials for the symmetric 
case, an optimal perturbation and the critical value, see Equation (j44p . 

In Figure [3] we present results for a three dimensional problem with the 
symmetric matrix 

5 = diag(l, 0.1, 0.01). (70) 

The spectral gap of the optimally perturbed nonreversible matrix (and of the 
generator of the semigroup) is given by 
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Fig. 3 Norms of the matrix exponentials for the 3x3 diagonal matrix (|70p and 
its optimal nonreversible perturbation. 
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t 



Fig. 4 Norms of the matrix exponentials for a diagonal matrix with random uni- 
formly distributed entries and its optimal nonreversible perturbation for A'^ — 100. 

which is a substantial improvement over that of S, namely 0.01. 

In Figurcdiwc consider a 100 x 100 diagonal matrix with random entries, 
uniformly distributed on [0, 1]. For our example the minimum diagonal ele- 
ment (spectral gap) is 0.0012. On the contrary, the spectral gap of Bj with 
J = Jopt is 0.4762. 

Finally, in Figure [5] we consider a drift that is a (high dimensional) fi- 
nite difference approximation of of the Laplacian with periodic boundary 
conditions. More precisely, consider the drift matrix 

Bii = 2, -Bi^i+i = -Bj-i,j = —1, 

with N ^ 100. In this case the improvement on the convergence rate is over 
three orders of magnitude, since 

TrS" 

min(cr(S)) = 9.67 x 10"'', whereas Re{(j{Bj)) = — = 2. 

Since the computational cost of calculating the optimal nonreversible per- 
turbation is very low, we believe that the algorithm developed in this paper 
can be used for sampling Gaussian distributions in infinite dimensions. 
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t 



Fig. 5 Norms of the matrix exponentials for the the discrete Laplacian and its 
optimal nonreversible perturbation for A'^ = 100. 

The algorithm developed in this paper provides us with the optimal non- 
reversible perturbation only in the case of linear drift. However, even for 
nonlinear problems it is always the case that the addition of a nonreversible 
perturbation can accelerate the convergence to equilibrium, as mentioned in 
the introduction. This is particularly the case for systems with metastable 
states and/or multiscale structure [25^; for such systems, a "clever" choice of 
the nonreversible perturbation can lead to a very significant increase in the 
rate of convergence to equilibrium. A systematic methodology for obtaining 
the optimal nonreversible perturbation for general reversible diffusions (i.e. 
not necessarily with a linear drift) will be developed elsewhere. 

We illustrate the advantage of adding a nonreversible perturbation to the 
dynamics by considering a few simple two-dimensional examples. In partic- 
ular, we consider the nonreversible dynamics 

dXt = {-I + SJ)VV{Xt)dt+^/2^dWt, (71) 

with (5 G M and J the standard 2x2 antisymmetric matrix, i.e. J12 = 
1, J21 = —1. For this class of nonreversible perturbations the parameter that 
we wish to choose in an optimal way is S. From our numerical experiments, 
we observed that even a non-optimal choice of S significantly accelerates con- 
vergence to equilibrium. To illustrate the effect of adding a nonreversible 
perturbation, we solve numerically (j7ip using the Euler-Marayama method 
with a sufficiently small time step and for a sufficiently large number of real- 
izations of the noise. We then compute the expectation value of observables 



42 



of the solution, in particular, the second moment by averaging over all the 
trajectories that we have generated. 

We use one of the potentials that were considered in [31], namely 

V{x,v) = \{x^~lf + \v\ (72) 



--■8 = 




D 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 



t 

Fig. 6 Second moment as a function of time for (jZip with the potential \T2} . We 
take as an initial condition and f3~^ — 0.1. 

In Figure [6] we present the convergence of the second moment to its equi- 
librium value for =0.1. Even in this very simple example, the addition 
of a nonreversible perturbation, with (5 = 10, speeds up convergence to equi- 
librium. Notice also that, as expected, the nonreversible perturbation leads 
to an oscillatory transient behavior. 

A Wick calculus 

In this article, we use a specific positivity property of the Wick calculus, which 
must not be confused with the more general and robust positivity property of the 
anti-Wick calculus This appendix recalls the basic facts about Wick calculus and 
its positivity property in L'^(R^, da;; C) . We refer the reader for details to 0(1] 
and references therein. This calculus is modelled on the creation and annihilation 
operators, a*j = —dxj + \xj and aj = dxj + of the A'^-dimensional harmonic 

oscillator Hamiltonian —Ax+^-^ -y ~ ^j'^j • The kernel of this Hamiltonian 

(the so-called vacuum state) is the Gaussian function 
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Let (ei, . . . , ejv) be an orthonormal basis of and use the notation 

The canonical commutation relations are (for two operators A and B, we denote 
[A, B] — AB — BA the commutator bracket) 

[a{ei),a{ej)] = [a* (e^), a* (ej)] =0, [a(e,), a' (e^)] ^ 5^J 
and become by sesquilinearity of the complex scalar product 

Vzi, 22 G , [a{zi),a*{z2)] = {zi , Z2)c , 
after setting for z £ C^, z = X^^i Z"'ej 

JV JV 

We use the multi-index notation introduced in Section 15.11 The orthonormal basis 
of Hermite functions in L'^(R^, da;; C) is (0„)„gj>iN given by 

=a*(ei)"i ...a*(ejv)""l«) = -^{a{e)T\Q). 



ni! . . . njv! Vw 



Since the basis of the subspace C(/>,i is by construction indexed by the n's in 

such that |n| — p , it can be considered as the p-fold symmetric tensor product 
of C^: 



with 5p(zi ® ■ • • ® 2p) = ^ ^ 2^(1 



The symmetrization operator Sp is actually the orthogonal projection from (C^)®^ 
onto V . This provides the description of L^(K^, dx\ C) as the so-called bosonic 



Fock space over 



oc p 

JV 



p=0 

where the infinite sum is orthogonal and complete. It is convenient to introduce 
also the algebraic orthogonal sum 

alg P l_^|2 

^ = V 'C'^ = 'CI*!' ■ • ■ ' ^iv]e-^ . 

pSN 

We now consider polynomials of 2N real coordinates {x, y) with z = x~>riy , written 
in the complex notation (z,z) as elements of 

c[F,...,i^,^\...,2^]= c,,pF,...,P^,z\...,^^], 

(p,q)eN2 
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where Cq,p[z^, . . . , z^, z^, . . . , 2^] denote the set of monomials homogeneous with 
degree q G N with respect to z and homogeneous with degree p G N with respect 
to z. Notice that monomials b G Cq,p[z^, . . . , z^ , z^ , . . . , z^] can be written 

b{z) = {z^" , te®^)c with 6 e L C"^; V C^^ . 

This provides a bijection between the sets Cq,p[z^, . . . . . . , z^] and L(V^ C^; V'' C^) 

with the inversion formula 

q'.pl 

With any monomial b G Cq.p[z^ , . . . , z^ , z''' , . . . , z^] (and by linearity with any 

polynomial of {z,'z)) we can associate an operator : V ^ V G {TSJ^ , dx;C) 

called its Wick quantization: 

When 6 G Cq,p [F, . . . , 2\ . . . , 2^] , and for any n G N its restriction b^'"'' | y„ : 
V" ^ \/"+q-p f^n+q-p jofincd by 



\yn cN = i-\p: + oo)[n) {n-py. ^ri+q-pib iSl lyn-p cN )i>„ ■ (76) 

Here are a few examples 

- if biz) = , 2)c , C G , then 6^-'= = a{0 ; 

- if b{z) = (z , C)c, C G , then b^"''' = o*(C) ; 

- if b{z) = (« , Az)c with yl G L{C^- C^) , then one recovers the second quantized 
version of A 

N 

tWick 



n-l 

— if 6(2;) = 0^=1 (-2^ 5 Cfe)c X nj=i(G > -*)c , so that the associated linear function 
infeG L(V''C^;V'C~) is 

b = <S,(|Ci) ® • • • (8) |C<i» ® <5p((6l ® • • • ® {(.p\) , 

and then 

fe^-'==a*(Ci)...a*(C,)a(ei)...a(ep). (74) 

For any polynomial 6 G C[z^ , . . . , z^ , z'^ , . . . , z'^] , fc G N and any z G C"'^, the 
fc-th order differential {d';b){z) is a C-lincar form on V'' while {d^b){z) is a 
C-antilinear form which can be identified with a vector via the complex scalar 
product 

k 

d^biz):u€\/C^ ^{u,d^b{z))c. 

We use the notation i.v for the C-bilinear duality product between £ G (V*' C")**^ 
and « G V*" . For any &i, &2 G C[F, . . . ,z^ , , . . . , z^], fc G N and all z G , 
the quantity 9^61(2) .(^62(2) is well defined in C and this defines a new polynomial 

d%.d^b2ec[^,...,7^,z\...,z'']. 
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Proposition 13 1. For any monomial b £ Cp^p[z^, . . . , z'^ , , . . . , z^\ such that 
6 > 0, the Wick quantized operator is non negative on V: 

V^GO, (<^,fo^^^V>L2 >0. 

2. When b e C[z^ , . . . ,7^ , z\ . . . , z'^] the formal adjomt of b''^'"'^ defined on V by 

\/ip,%l) £V, {if, {b^"'^)'^)L^ = (b^'^^ip , '!/')l2 

IS given by {b^^'^'')' = (6) . 

3. The set of polynomials C[z-'-, . . . , z^ , z^ , . . . , z'^] is an algebra for the operation 

f i_ iiWickj \Wick j Wick jWick /t-» , /r-* 
(,OlP 02) ~ bi O O2 : V ^ V 

With biC^%2 = J2i^9%.d^b2, (75) 

fc=0 

where the sum in the right-hand side is actually finite. 

Proof 1) It comes from the definition H73|l with p = q . Actually for any n G N, 
b (g) lyn f^N is non negative and 

Sn (b ® Jy^ 5„ = 5,t (b ® /yn c") <S„ > . 

For = ©J^"q^ vp„ G D we get 

n — 

2) It results from the definition (I73|l after noticing 

6(1) = {z®i , bz®p)c = (z®" , 6*z®'')c, 
where 6* is the adjoint of 6 G -C(V'' C^; V*" C^) . 

3) The definition (|73|) ensures that for & G C[z^ , . . . , 2^, 2:^, . . . , z^] , the operator 
b^^"'' sends D into itself, so that bf'"^ ob'^^"'' is well defined. By linearity the result 
comes from considering the specific case when bf'"^ and bf^"^ have the form ([7^1). 
With the polarization identities 



5,(|Ci)®---®iC.)) = 
5j,({Cii®---«){ep|) = 




the problem is reduced to bf""'' = [a* (C)]"' [a(r)]''' for i = 1,2. But this is a 
simple iterated application of [a(5^), a*((^^)] — (^^ , C^)c- □ 
A useful consequence of the properties stated above for our case is the following 
lemma. It provides a lower bound for a differential operator with a specific quartic 
symbol. Of course, it is in a very specific case but it is much stronger than what 
would give the Feffermann-Phong inequality. The refore, it is probably not easily 
accessible via the Weyl or anti-Wick calculus (see [221127) '). 
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Lemma 4 Let S,Q, J be real matrices such that S £ 5^''(K), Q G 5^''(R) and J G 
.4jv(K). Let us consider the operator L = ~Cj (see Equation \2T\ for the definition) 
and C = Cq (see Equation (I60p j. The operator L ( respectively C) is the Wick 
quantization of the polynomial l(z) = (z , Sz)c — (z , Jz)c (respectively p{z) = 
(z , Qz)c). Moreover, we have the following estimate: \/ip £ V, 

(¥3, (L*C + CL)(^)i2 > (|^, ((z, [50 + QS'+ JQ-Qj] z)^)'^"%^ . (76) 

Proof The fact that L = (^(z))^*=* and C = {p{z))^''''' is easy to check. Notice 
that the polynomials £ and p satisfy 

Re £(z) = (z , 5z)c , Im /(z) = — i(z , Jz)c , p(z) = p(z) . 

I 

We are looking for a lower bound for L*C + CL. Using formula (|75p . the Wick 
symbol of L'C + CL is 

I(z)p{z) + p(z) ^(z) + a.I(iy . chp(z) + 9,p(z) . d^l{z) 

= (z®% (5'®0 + 0®S)z®2)^+ (z, (S-Q + Q5 + JQ - j) z)^. 

Since S and Q are non negative matrices, we deduce that S®Q and Q® S are non 
negative an d th e first term is thu s non negative. By applying the first statement of 
Proposition 1131 one obtains (|76p . □ 
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